One-dimensional finite element method
For a given computational mesh \( \{ T_{hp }= \left(K, X\left(K\right), \Pi_p \right) \}_K \) a reference mesh is a family \( T_{\frac{h}{2}p+1} \) of finite elements
\( T_{\frac{h}{2}p+1}=\{ \left( K, X\left(K\right), \Pi_p \right) \}_K { } \) such that
\( \forall K \in T_{ \frac{h}{2}p+1} \exists K_1,K_2 \in {\cal P}(T_{hp },K) \) such that
\( K=K_1 \cup K_2, \textrm{meas}K_1\cap K_2=0 \),
\( X(K_1),X(K_2) \in {\cal P}(T_{hp },X(K)), dim X(K)=dimX(K_1)+1=dimX(K_2)+1 \) where \( {\cal P } (T_{hp},K ) \) and \( {\cal P }(T_{hp},X(K)) \) denotes projections onto the first and the second component \( \left( K, X\left( K \right), \Pi_p \right) \).
Approximation space over a one-dimensional computational mesh \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) is defined as
\( V_{hp }= \{ v \in C(\Omega): \forall K \in {\cal P }( T_{hp},K): {\cal P }( v,K) \in X(K) \} \)
where \( {\cal P }( T_{hp},K) \) is a set of intervals representing the geometry of elements drawn from a triplet representing a one-dimensional computational mesh, \( {\cal P }( v,K) \) is the projection of a function onto a compartment representing the geometry of an element.
Let \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) denote one-dimensional computational mesh.
Let \( \{ e_i^{hp} \}_i \) denote the basis \( V_{hp }= span\{ e_i^{hp} \} \).
Let \( \chi^K_k \in X \left( K\right) \) denote the shape functions over element \( K \).
Then \( \forall K \in {\cal P }( T_{hp},K), \forall i, \exists k : {\cal P}(e_i^{hp},K) = \chi^K_k \).
there exists an inverse operator \( {\cal I}^2 \ni (k,K)\rightarrow i(k,K)\in {\cal I} \) that rewrites
\( i(k,K) \) the index of i-th global basis function related to \( k \)-th local shape function over an element \( K \).
In the cases considered in this manual, this mapping is isomorphic.
Approximation space over one-dimensional reference mesh \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) is defined as
\( V_{\frac{h}{2}p+1 }= \{ v \in C(\Omega): \forall K \in {\cal P }( T_{\frac{h}{2}p+1},K): {\cal P }( v,K) \in X(K) \} \)
where \( {\cal P }( T_{\frac{h}{2}p+1},K) \) is a set of intervals representing the geometry of elements drawn from a triple representing a one-dimensional computational mesh, \( {\cal P }( v,K) \) is the projection of a function onto a compartment representing the geometry of an element.
Let \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) denote one-dimensional reference mesh.
Let \( \{ e_i^{\frac{h}{2}p+1} \}_i \) denote the basis of space \( V_{\frac{h}{2}p+1 }= span\{ e_i^{\frac{h}{2}p+1} \} \).
Let \( \chi^K_k \in X \left( K\right) \) denote shape functions over an element \( K \).
Then \( \forall K \in {\cal P }( T_{\frac{h}{2}p+1},K), \forall i, \exists k : {\cal P}(e_i^{\frac{h}{2}p+1},K) = \chi^K_k \).
There exists and inverse operator \( {\cal I}^2 \ni (k,K)\rightarrow i(k,K)\in {\cal I} \) that assigns
\( i(k,K) \) the index of global i-th basis functions related to \( k \)-th local shape function over an element \( K \).
In the cases considered in this manual, this mapping is isomorphic.
The reference mesh is used to estimate the relative error on the computational mesh. Often, when a computational mesh is described in the context of a reference mesh, it is called a coarse mesh, and the reference mesh is called a fine mesh, because it is created by breaking elements of the coarse mesh into smaller elements and increasing the degree of polynomial approximation by one.
For a given computational mesh \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) find coefficients \( \{ u^{hp}_i \}_{i=1,...,N^{hp}} \) of an approximate solution \( V \supset V_{hp } \ni u_{hp }=\sum_{i=1,...,N^{hp } } u_i^{hp } e_i^{hp } \) such that \( \sum_{i=m,...,N^{hp } } u_m^{hp } B(e_m^{hp},e_n^{hp})=L(e_n^{hp}), n=1,...,N^{hp} \) where \( B(e_m^{hp},e_n^{hp})=\int_0^l \left( a(x) \frac{de_m^{hp}(x)}{dx} \frac{de_n^{hp}(x)}{dx} +b(x)\frac{de_m^{hp}(x)}{dx}e_n^{hp}(x)+c(x)e_m^{hp}(x)e_n^{hp}(x)\right)dx + \beta e_m^{hp}(l)e_n^{hp}(l) \) and \( L(e_n^{hp})=\int_0^l f(x) e_n^{hp}(x)dx + \gamma e_n^{hp}(l) \). The basis \( \{ e^{hp}_i \}_{i=1,...,N^{hp}} \) of approximation space \( V_{hp} \) is obtained by merging (glueing together) shape functions from space \( X\left(K\right) = {\cal P}(T_{hp},X(K)) \) for particular finite elements from the computational mesh \( T_{hp} \) into the global basis functions.
For a given computational mesh \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) find coefficients \( \{ u^{\frac{h}{2}p+1}_i \}_{i=1,...,n^{\frac{h}{2}p+1}} \) of an approximate solution \( V \supset V_{\frac{h}{2}p+1 } \ni u_{\frac{h}{2}p+1 }=\sum_{i=1,...,N^{\frac{h}{2}p+1 } } u_i^{\frac{h}{2}p+1 } e_i^{\frac{h}{2}p+1 } \) such that \( \sum_{m=1,...,N^{\frac{h}{2}p+1 } } u_m^{\frac{h}{2}p+1 } B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1})=L(e_n^{\frac{h}{2}p+1}), n=1,...,N^{hp} \) where \( B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1})=\int_0^l \left( a(x) \frac{de_m^{\frac{h}{2}p+1}(x)}{dx} \frac{de_n^{\frac{h}{2}p+1}(x)}{dx} +b(x)\frac{de_m^{\frac{h}{2}p+1}(x)}{dx}e_n^{\frac{h}{2}p+1}(x)+c(x)e_m^{\frac{h}{2}p+1}(x)e_n^{\frac{h}{2}p+1}(x)\right)dx +\\+ \beta e_m^{\frac{h}{2}p+1}(l)e_n^{hp}(l) \) and \( L(e_n^{\frac{h}{2}p+1})=\int_0^l f(x) e_n^{\frac{h}{2}p+1}(x)dx + \gamma e_n^{\frac{h}{2}p+1}(l) \).
The basis \( \{ e^{\frac{h}{2}p+1}_i \}_{i=1,...,N^{\frac{h}{2}p+1}} \) of an approximation space \( V_{hp} \) is obtained by merging (glueing together) shape functions from space \( X\left(K\right) = {\cal P}(T_{\frac{h}{2}p+1},X(K)) \) for particular elements from the computational mesh \( T_{\frac{h}{2}p+1} \) into the global basis functions.
Coefficients \( \{ u^{hp}_i \}_{i=1,...,N^{hp}} \) are calculated by solving a system of equations
\( {\bf A}{\bf u}={\bf b} \) where
\( {\bf A}= \begin{bmatrix} B(e_1^{hp},e_1^{hp}) & \cdots & B(e_{N^{hp}}^{hp},e_1^{hp}) \\ \vdots & \cdots & \vdots \\ B(e_1^{hp},e_{N^{hp}}^{hp}) & \cdots & B(e_{N^{hp}}^{hp},e_{N^{hp}}^{hp}) \end{bmatrix} \).
\( {\bf u} = \begin{bmatrix} u_1^{hp} \\ \vdots\\ u_{N^{hp}}^{hp} \end{bmatrix} \).
\( {\bf b} = \begin{bmatrix} L(e_1^{hp}) \\ \vdots \\ L(e_{N^{hp}}^{hp}) \end{bmatrix} \).
The Dirichlet boundary condition is obtained by resetting the rows corresponding to the basis functions defined on the Dirichlet boundary, placing the one on the diagonal, and resetting the right-hand side. The effect of the boundary condition on our system of equations was transferred to the right-hand side by the shift operation.
Integrals \( B(e_m^{hp},e_n^{hp}) \) and \( L(e_n^{hp}) \) are calculated element by element by using Gaussian quadratures.
\( B(e_m^{hp},e_n^{hp})=\sum_k \left( a(x^k) \frac{de_m^{hp}(x^k)}{dx} \frac{de_n^{hp}(x^k)}{dx} +b(x^k)\frac{de_m^{hp}(x^k)}{dx}e_n^{hp}(x^k)+c(x^k)e_m^{hp}(x^k)e_n^{hp}(x^k)\right)Jac(x^k)w^k +\\+ \beta e_m^{hp}(l)e_n^{hp}(l) \) and \( L(e_n^{hp})=\sum_k f(x^k) e_n^{hp}(x^k) Jac(x^k)w^k + \gamma e_n^{hp}(l) \)
where \( Jac(x^k) \) is the value at the Gaussian quadrature point of the Jacobian mapping from the master element to the element.
Generation of the system of equations for a given computational mesh \( T_{hp}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \) can be obtained by using the following algorithm:
1 for \( K \in {\cal P}(T_{hp},K) \) (loop through finite element mesh)
2 for \( \psi_k \in X(K) \) (loop through shape functions from element \( K \) )
3 \( {\bf b}(i(k,K)) += L(\psi_k) \) (calculation of \( K \) element contribution to the integral)
4 for \( \psi_l \in X(K) \) (loop through shape functions from element \( K \) )
5 \( {\bf A}(i(k,K),i(l,K)) += B(\psi_k,\psi_l) \) (calcuiation of \( K \) element contribution to the integral).
Coefficients \( \{ u^{\frac{h}{2}p+1}_i \}_{i=1,...,n^{\frac{h}{2}p+1}} \) are calculated by solving a system of equations \( {\bf A }{\bf u} = {\bf b} \)
\( {\bf A} = \begin{bmatrix} B(e_1^{\frac{h}{2}p+1},e_1^{\frac{h}{2}p+1}) & \cdots & B(e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1},e_1^{\frac{h}{2}p+1}) \\ \vdots & \cdots & \vdots \\ B(e_1^{\frac{h}{2}p+1},e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1}) & \cdots & B(e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1},e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1}) \end{bmatrix} \).
\( {\bf u} = \begin{bmatrix} u_1^{\frac{h}{2}p+1} \\ \vdots\\ u_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1} \end{bmatrix} \).
\( {\bf b} =\begin{bmatrix} L(e_1^{\frac{h}{2}p+1}) \\ \vdots \\ L(e_{N^{\frac{h}{2}p+1}}^{\frac{h}{2}p+1}) \end{bmatrix} \).
Integrals \( B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1}) \) and \( L(e_n^{hp}) \) are computed element by element using the Gaussian quadratures.
\( B(e_m^{\frac{h}{2}p+1},e_n^{\frac{h}{2}p+1})=\sum_k \left( a(x^k) \frac{de_m^{\frac{h}{2}p+1}(x^k)}{dx} \frac{de_n^{\frac{h}{2}p+1}(x^k)}{dx} +b(x^k)\frac{de_m^{\frac{h}{2}p+1}(x^k)}{dx}e_n^{\frac{h}{2}p+1}(x^k)+c(x^k)e_m^{\frac{h}{2}p+1}(x^k)e_n^{\frac{h}{2}p+1}(x^k)\right) *\\ *Jac(x^k)w^k + \beta e_m^{\frac{h}{2}p+1}(l)e_n^{\frac{h}{2}p+1}(l) \) and \( L(e_n^{\frac{h}{2}p+1})=\sum_k f(x^k) e_n^{\frac{h}{2}p+1}(x^k) Jac(x^k)w^k + \gamma e_n^{\frac{h}{2}p+1}(l) \)
where \( Jac(x^k) \) is the value at the Gaussian quadrature point of the Jacobian mapping from the master element to the element.
Generation of the system of equations for a given computational mesh \( T_{\frac{h}{2}p+1}=\{ \left(K, X\left(K\right), \Pi_p \right) \}_K \)
can be obtained by using the following algorithm:
1 for \( K \in {\cal P}(T_{\frac{h}{2}p+1 },K) \) (loop through elements of the computational meshj)
2 for \( \psi_k \in X(K) \) (loop through shape functions of element \( K \) )
3 \( {\bf b}(i(k,K)) += L(\psi_k) \) (loop through \( K \) element contribution to the integral)
4 for \( \psi_l \in X(K) \) (loop through shape functions over an element \( K \) )
5 \( {\bf A}(i(k,K),i(l,K)) += B(\psi_k,\psi_l) \) (computations of \( K \) element contribution to the integral).